R Notebook

library(Seurat)
library(scHolography)
library(dplyr)
scHolography.obj <- readRDS("~/OneDrive - Northwestern University/Northwestern/Yi Lab/Manuscript/Code_2ndDraft/Human Skin/scHolography.obj_new.rds")

sup.with.more.basal <- which(rowSums(scHolography.obj$adj.mtx[which(scHolography.obj$scHolography.sc$celltype%in%"Suprabasal"),which(scHolography.obj$scHolography.sc$celltype%in%"Basal")])>=min(boxplot.stats(rowSums(scHolography.obj$adj.mtx[which(scHolography.obj$scHolography.sc$celltype%in%"Suprabasal"),which(scHolography.obj$scHolography.sc$celltype%in%"Basal")]))$out)
)

sup.with.more.basal.names <- (colnames(scHolography.obj$scHolography.sc)[which(scHolography.obj$scHolography.sc$celltype%in%"Suprabasal")])[sup.with.more.basal]
sup.with.less.basal.names <- (colnames(scHolography.obj$scHolography.sc)[which(scHolography.obj$scHolography.sc$celltype%in%"Suprabasal")])[-sup.with.more.basal]

scHolography.obj$scHolography.sc$suprabasal.sub <- as.character(scHolography.obj$scHolography.sc$celltype)
scHolography.obj$scHolography.sc$suprabasal.sub[sup.with.more.basal.names] <- "Transition_KC"
scHolography.obj$scHolography.sc$suprabasal.sub[sup.with.less.basal.names] <- "Differentiated_KC"
scHolography.obj$scHolography.sc$suprabasal.sub <-factor(as.character(scHolography.obj$scHolography.sc$suprabasal.sub ),levels = c("Basal", "Transition_KC", "Differentiated_KC","Endothelial", "Fibroblast","Smooth Muscle", "Lymphatic Endothelial","Schwann","Glandular Epithelium","Immune","Melanocyte"))

scHolography.obj$scHolography.sc <- SetIdent(scHolography.obj$scHolography.sc,value = "suprabasal.sub") 
DefaultAssay(scHolography.obj$scHolography.sc) <- "RNA"
scHolography.obj$scHolography.sc <- ScaleData(scHolography.obj$scHolography.sc)

Centering and scaling data matrix

  |                                                                                                                                                
  |                                                                                                                                          |   0%
  |                                                                                                                                                
  |======                                                                                                                                    |   4%
  |                                                                                                                                                
  |============                                                                                                                              |   9%
  |                                                                                                                                                
  |==================                                                                                                                        |  13%
  |                                                                                                                                                
  |========================                                                                                                                  |  17%
  |                                                                                                                                                
  |==============================                                                                                                            |  22%
  |                                                                                                                                                
  |====================================                                                                                                      |  26%
  |                                                                                                                                                
  |==========================================                                                                                                |  30%
  |                                                                                                                                                
  |================================================                                                                                          |  35%
  |                                                                                                                                                
  |======================================================                                                                                    |  39%
  |                                                                                                                                                
  |============================================================                                                                              |  43%
  |                                                                                                                                                
  |==================================================================                                                                        |  48%
  |                                                                                                                                                
  |========================================================================                                                                  |  52%
  |                                                                                                                                                
  |==============================================================================                                                            |  57%
  |                                                                                                                                                
  |====================================================================================                                                      |  61%
  |                                                                                                                                                
  |==========================================================================================                                                |  65%
  |                                                                                                                                                
  |================================================================================================                                          |  70%
  |                                                                                                                                                
  |======================================================================================================                                    |  74%
  |                                                                                                                                                
  |============================================================================================================                              |  78%
  |                                                                                                                                                
  |==================================================================================================================                        |  83%
  |                                                                                                                                                
  |========================================================================================================================                  |  87%
  |                                                                                                                                                
  |==============================================================================================================================            |  91%
  |                                                                                                                                                
  |====================================================================================================================================      |  96%
  |                                                                                                                                                
  |==========================================================================================================================================| 100%

scHolographyPlot(scHolography.obj,color.by = "suprabasal.sub",palette = "Set3")

my.col <- c("#A6CEE3","#B96A59","#743A36")
scHolography.obj.sub<- scHolography.obj
scHolography.obj.sub$scHolography.sc <- subset(scHolography.obj.sub$scHolography.sc,cells=which(scHolography.obj.sub$scHolography.sc$suprabasal.sub%in%c("Basal", "Transition_KC", "Differentiated_KC")))
scHolography.obj.sub$scHolography.sc$suprabasal.sub <- droplevels(scHolography.obj.sub$scHolography.sc$suprabasal.sub)


scene = list(camera = list(eye = list(x = -1, z = 0, y = 2)))
fig1 <- scHolographyPlot(scHolography.obj.sub,color.by = "suprabasal.sub",color = my.col)%>% plotly::layout(scene = scene)
fig1

tiff("~/OneDrive - Northwestern University/Northwestern/Yi Lab/Manuscript/Code_3rdDraft//sup.num.basal.vs.sup.neighbor.tiff", units="in", width=5, height=5, res=300)
ggplot(data = sup.epi.neighb, aes(x = n.suprabasal, y = n.basal)) +
  geom_jitter(aes(color=type)) +
  stat_cor(label.y = 30, size = 6) +theme_classic()+scale_color_manual(values = my.col[2:3])+ 
  geom_smooth(method = "lm", 
              formula = y ~ x) +NoLegend()+xlab("# of Suprabasal Neighbors")+ylab("# of Basal Neighbors")+theme(axis.text = element_text(size = 15,face = "bold"),axis.text.x = element_text(angle = 45,hjust=1),axis.title=element_text(size = 15,face = "bold"))
dev.off()

null device 
          1 

my_comparisons <- list(  c("Basal", "Transition_KC"), c("Transition_KC", "Differentiated_KC")  )

#tiff("~/OneDrive - Northwestern University/Northwestern/Yi Lab/Manuscript/Code_2ndDraft/Human Skin//clusterDist.Basal.tiff", units="in", width=5, height=5, res=300)
scHolography::clusterDistanceBoxplot(scHolography.obj.sub,annotationToUse = "suprabasal.sub",query.cluster.list = c("Basal", "Transition_KC","Differentiated_KC"),reference.cluster = c("Basal"))+ggplot2::scale_fill_manual(values =my.col)+ 
  geom_signif(comparisons = my_comparisons,test.args = list(size=5),map_signif_level = TRUE)+NoLegend()+theme(axis.text = element_text(size = 16,face = "bold"),axis.title.x = element_blank(),axis.title.y = element_text(size = 20,face = "bold") ,axis.text.x = element_text(angle = 45,hjust=1))
#dev.off()

mark.ba <- FindMarkers(scHolography.obj$scHolography.sc,ident.1 =  "Transition_KC",ident.2 = "Basal" ,assay = "RNA")

  |                                                  | 0 % ~calculating  
  |+                                                 | 1 % ~01s          
  |++                                                | 2 % ~01s          
  |++                                                | 3 % ~01s          
  |+++                                               | 4 % ~01s          
  |+++                                               | 5 % ~01s          
  |++++                                              | 6 % ~01s          
  |++++                                              | 7 % ~01s          
  |+++++                                             | 8 % ~01s          
  |+++++                                             | 9 % ~01s          
  |++++++                                            | 10% ~01s          
  |++++++                                            | 11% ~01s          
  |+++++++                                           | 12% ~01s          
  |+++++++                                           | 13% ~01s          
  |++++++++                                          | 14% ~01s          
  |++++++++                                          | 15% ~01s          
  |+++++++++                                         | 16% ~01s          
  |+++++++++                                         | 17% ~01s          
  |++++++++++                                        | 18% ~01s          
  |++++++++++                                        | 19% ~01s          
  |+++++++++++                                       | 20% ~01s          
  |+++++++++++                                       | 21% ~01s          
  |++++++++++++                                      | 22% ~01s          
  |++++++++++++                                      | 23% ~01s          
  |+++++++++++++                                     | 24% ~01s          
  |+++++++++++++                                     | 26% ~01s          
  |++++++++++++++                                    | 27% ~01s          
  |++++++++++++++                                    | 28% ~01s          
  |+++++++++++++++                                   | 29% ~01s          
  |+++++++++++++++                                   | 30% ~01s          
  |++++++++++++++++                                  | 31% ~01s          
  |++++++++++++++++                                  | 32% ~01s          
  |+++++++++++++++++                                 | 33% ~01s          
  |+++++++++++++++++                                 | 34% ~01s          
  |++++++++++++++++++                                | 35% ~01s          
  |++++++++++++++++++                                | 36% ~01s          
  |+++++++++++++++++++                               | 37% ~01s          
  |+++++++++++++++++++                               | 38% ~01s          
  |++++++++++++++++++++                              | 39% ~01s          
  |++++++++++++++++++++                              | 40% ~01s          
  |+++++++++++++++++++++                             | 41% ~01s          
  |+++++++++++++++++++++                             | 42% ~01s          
  |++++++++++++++++++++++                            | 43% ~01s          
  |++++++++++++++++++++++                            | 44% ~01s          
  |+++++++++++++++++++++++                           | 45% ~01s          
  |+++++++++++++++++++++++                           | 46% ~01s          
  |++++++++++++++++++++++++                          | 47% ~01s          
  |++++++++++++++++++++++++                          | 48% ~01s          
  |+++++++++++++++++++++++++                         | 49% ~01s          
  |+++++++++++++++++++++++++                         | 50% ~01s          
  |++++++++++++++++++++++++++                        | 51% ~01s          
  |+++++++++++++++++++++++++++                       | 52% ~00s          
  |+++++++++++++++++++++++++++                       | 53% ~00s          
  |++++++++++++++++++++++++++++                      | 54% ~00s          
  |++++++++++++++++++++++++++++                      | 55% ~00s          
  |+++++++++++++++++++++++++++++                     | 56% ~00s          
  |+++++++++++++++++++++++++++++                     | 57% ~00s          
  |++++++++++++++++++++++++++++++                    | 58% ~00s          
  |++++++++++++++++++++++++++++++                    | 59% ~00s          
  |+++++++++++++++++++++++++++++++                   | 60% ~00s          
  |+++++++++++++++++++++++++++++++                   | 61% ~00s          
  |++++++++++++++++++++++++++++++++                  | 62% ~00s          
  |++++++++++++++++++++++++++++++++                  | 63% ~00s          
  |+++++++++++++++++++++++++++++++++                 | 64% ~00s          
  |+++++++++++++++++++++++++++++++++                 | 65% ~00s          
  |++++++++++++++++++++++++++++++++++                | 66% ~00s          
  |++++++++++++++++++++++++++++++++++                | 67% ~00s          
  |+++++++++++++++++++++++++++++++++++               | 68% ~00s          
  |+++++++++++++++++++++++++++++++++++               | 69% ~00s          
  |++++++++++++++++++++++++++++++++++++              | 70% ~00s          
  |++++++++++++++++++++++++++++++++++++              | 71% ~00s          
  |+++++++++++++++++++++++++++++++++++++             | 72% ~00s          
  |+++++++++++++++++++++++++++++++++++++             | 73% ~00s          
  |++++++++++++++++++++++++++++++++++++++            | 74% ~00s          
  |++++++++++++++++++++++++++++++++++++++            | 76% ~00s          
  |+++++++++++++++++++++++++++++++++++++++           | 77% ~00s          
  |+++++++++++++++++++++++++++++++++++++++           | 78% ~00s          
  |++++++++++++++++++++++++++++++++++++++++          | 79% ~00s          
  |++++++++++++++++++++++++++++++++++++++++          | 80% ~00s          
  |+++++++++++++++++++++++++++++++++++++++++         | 81% ~00s          
  |+++++++++++++++++++++++++++++++++++++++++         | 82% ~00s          
  |++++++++++++++++++++++++++++++++++++++++++        | 83% ~00s          
  |++++++++++++++++++++++++++++++++++++++++++        | 84% ~00s          
  |+++++++++++++++++++++++++++++++++++++++++++       | 85% ~00s          
  |+++++++++++++++++++++++++++++++++++++++++++       | 86% ~00s          
  |++++++++++++++++++++++++++++++++++++++++++++      | 87% ~00s          
  |++++++++++++++++++++++++++++++++++++++++++++      | 88% ~00s          
  |+++++++++++++++++++++++++++++++++++++++++++++     | 89% ~00s          
  |+++++++++++++++++++++++++++++++++++++++++++++     | 90% ~00s          
  |++++++++++++++++++++++++++++++++++++++++++++++    | 91% ~00s          
  |++++++++++++++++++++++++++++++++++++++++++++++    | 92% ~00s          
  |+++++++++++++++++++++++++++++++++++++++++++++++   | 93% ~00s          
  |+++++++++++++++++++++++++++++++++++++++++++++++   | 94% ~00s          
  |++++++++++++++++++++++++++++++++++++++++++++++++  | 95% ~00s          
  |++++++++++++++++++++++++++++++++++++++++++++++++  | 96% ~00s          
  |+++++++++++++++++++++++++++++++++++++++++++++++++ | 97% ~00s          
  |+++++++++++++++++++++++++++++++++++++++++++++++++ | 98% ~00s          
  |++++++++++++++++++++++++++++++++++++++++++++++++++| 99% ~00s          
  |++++++++++++++++++++++++++++++++++++++++++++++++++| 100% elapsed=01s  

# VlnPlot(scHolography.obj$scHolography.sc,c("ITGB4","LAMB3","DST","COL17A1"),idents = c("Transition_KC","Differentiated_KC","Basal "),assay = "SCT",ncol = 2,cols = my.col)
# VlnPlot(scHolography.obj$scHolography.sc,c("FABP7","THEM5","KRT10","DSG1"),idents = c("Transition_KC","Differentiated_KC","Basal"),assay = "SCT",ncol = 2,cols = my.col)

VlnPlot(scHolography.obj$scHolography.sc,c("KRT5"),idents = c("Transition_KC","Differentiated_KC","Basal"),assay = "SCT", cols = my.col)+ 
  geom_signif(comparisons = my_comparisons,test.args = list(size=5),map_signif_level = TRUE)+ylim(c(0,7))+NoLegend()+ theme(axis.title.x = element_blank(),axis.title.y = element_blank(),axis.text = element_text(size = 15,face = "bold"))

Scale for y is already present.
Adding another scale for y, which will replace the existing scale.

epi.sub <- subset(scHolography.obj$scHolography.sc, cells = which(scHolography.obj$scHolography.sc@active.ident%in% c("Transition_KC","Differentiated_KC")))
DefaultAssay(epi.sub) <- "RNA"
DotPlot(epi.sub,features = c("HSPG2","COL17A1","PDGFA","ITGB1","ITGA2","COL7A1","ITGB4","LAMB3","DST", 
                             "KRT1","KRT10","DSG1","EVPL","PERP", #keratinization
                             "PLA2G4F","STARD10","PLPP2","RORA","HMGCR","FABP7","DBI","THEM5","EPHX2","PTGS1","GSTA4" ),assay = "SCT")+viridis::scale_color_viridis()+ggplot2::coord_flip()+ theme(axis.text.x = element_text(angle = 45, hjust=1))

tiff("~/OneDrive - Northwestern University/Northwestern/Yi Lab/Manuscript/Code_2ndDraft/Human Skin//DotPlot.trans.vs.diff.tiff", units="in", width=3, height=10, res=300)
DotPlot(epi.sub,features = c(rownames(mark.sup %>%top_n(n = 10, wt = avg_log2FC) ),rownames(mark.sup %>%top_n(n = 10, wt = -avg_log2FC) )),assay = "SCT")+viridis::scale_color_viridis()+ggplot2::coord_flip()+ theme(axis.text.x = element_text(angle = 45, hjust=1))+theme(axis.title.x = element_blank(),axis.title.y = element_blank(),axis.text = element_text(size = 15,face = "bold"))+NoLegend()
dev.off()

epi.sub <- subset(scHolography.obj$scHolography.sc, cells = which(scHolography.obj$scHolography.sc@active.ident%in% c("Transition_KC","Basal")))
DefaultAssay(epi.sub) <- "RNA"
DotPlot(epi.sub,features = c("ANXA1","CLSTN1","DUSP6","PTPRZ1","MEF2A", #signal transduction
                             "COL17A1","PDGFA","ITGB1","COL4A5","COL7A1","ITGB4","LAMB3","DST", #extracellular matrix
                             "DSG1","KRT1","KRT10",
                             "NECTIN1","NECTIN4","CLDN4","CLDN1"
),assay = "SCT")+viridis::scale_color_viridis()+ggplot2::coord_flip()+ theme(axis.text.x = element_text(angle = 45, hjust=1))

tiff("~/OneDrive - Northwestern University/Northwestern/Yi Lab/Manuscript/Code_2ndDraft/Human Skin//DotPlot.trans.vs.basal.tiff", units="in", width=3, height=10, res=300)
DotPlot(epi.sub,features = c(rownames(mark.ba %>%top_n(n = 10, wt = avg_log2FC) ),rownames(mark.ba %>%top_n(n = 10, wt = -avg_log2FC) )),assay = "SCT")+viridis::scale_color_viridis()+ggplot2::coord_flip()+ theme(axis.text.x = element_text(angle = 45, hjust=1))+theme(axis.title.x = element_blank(),axis.title.y = element_blank(),axis.text = element_text(size = 15,face = "bold"))+NoLegend()
dev.off()

library(CellChat)
epi.sub <- subset(scHolography.obj$scHolography.sc, cells = which(scHolography.obj$scHolography.sc@active.ident%in% c("Transition_KC","Basal","Differentiated_KC")))
epi.sub$suprabasal.sub <- droplevels(epi.sub$suprabasal.sub )
cellchat.obj <- createCellChat(epi.sub,group.by = "suprabasal.sub",assay = "SCT")
CellChatDB <- CellChatDB.human
CellChatDB.use <- CellChatDB 
cellchat.obj@DB <- CellChatDB.use
cellchat.obj <- subsetData(cellchat.obj)
cellchat.obj <- identifyOverExpressedGenes(cellchat.obj)
cellchat.obj <- identifyOverExpressedInteractions(cellchat.obj)

cellchat.obj <- computeCommunProb(cellchat.obj)
# Filter out the cell-cell communication if there are only few number of cells in certain cell groups
cellchat.obj <- filterCommunication(cellchat.obj, min.cells = 10)
cellchat.obj <- computeCommunProbPathway(cellchat.obj)
cellchat.obj <- aggregateNet(cellchat.obj)

library(RColorBrewer)

cellchat.obj <- netAnalysis_computeCentrality(cellchat.obj, slot.name = "netP") 
my.sig <- c("LAMININ","THBS","DESMOSOME","NECTIN","NOTCH","EGF","NRG")
tiff("~/OneDrive - Northwestern University/Northwestern/Yi Lab/Manuscript/Code_2ndDraft/Human Skin//Epithelial.cellchat.out.tiff", units="in", width=9, height=6, res=300)
netAnalysis_signalingRole_heatmap(cellchat.obj, pattern = "outgoing",width = 9,height = 6,font.size = 15,color.heatmap = "YlOrRd",signaling =my.sig,color.use = my.col )
dev.off()

tiff("~/OneDrive - Northwestern University/Northwestern/Yi Lab/Manuscript/Code_2ndDraft/Human Skin//Epithelial.cellchat.in.tiff", units="in", width=9, height=6, res=300)
netAnalysis_signalingRole_heatmap(cellchat.obj, pattern = "incoming",width = 9,height = 6,font.size = 15,color.heatmap = "YlOrRd",signaling =my.sig,color.use = my.col  )
dev.off()

Reactome_2022_high_diff <- read.delim("~/Library/CloudStorage/OneDrive-NorthwesternUniversity/Northwestern/Yi Lab/Manuscript/Code_2ndDraft/Human Skin/Reactome_2022_HIGH_in_diff_kc.txt", quote="")
Reactome_2022_high_diff.sub <- Reactome_2022_high_diff[1:5,]
Reactome_2022_high_diff.sub$Term <- unlist(lapply(Reactome_2022_high_diff.sub$Term , function(x) strsplit(x," R-HSA")[[1]][1]))
Reactome_2022_high_diff.sub$Term <- factor(Reactome_2022_high_diff.sub$Term,levels = rev(Reactome_2022_high_diff.sub$Term))

tiff("~/OneDrive - Northwestern University/Northwestern/Yi Lab/Manuscript/Code_2ndDraft/Human Skin//Reactome_2022_high_diff.tiff", units="in", width=6, height=6, res=300)
ggplot(data=Reactome_2022_high_diff.sub, aes(x=Term, y=-log(Adjusted.P.value))) +
  geom_bar(stat="identity",fill=my.col[3]) + coord_flip()+theme_classic()+theme(axis.title.y = element_blank(),axis.text.y = element_blank(),axis.text = element_text(size = 15,face = "bold"),axis.title  = element_text(size = 15,face = "bold"))
dev.off()



ggplot(data=Reactome_2022_high_diff.sub, aes(x=Term, y=-log(Adjusted.P.value))) +
  geom_bar(stat="identity",fill=my.col[3]) + coord_flip()+theme_classic()+theme(axis.title.y = element_blank(),axis.text = element_text(size = 15,face = "bold"),axis.title  = element_text(size = 15,face = "bold"))+geom_text(aes(label=Genes), y=c(0),hjust = 0,color="white", size=c(1.75,2.25,2.5,2.5,2.5))+theme(axis.title.y = element_blank())

### Run Reatome 2022 on genes with <0.05 adjusted p values in UpReg_in_Differentiated_KC_comp_w_Transition_KC.csv and UpReg_in_Transition_KC_comp_w_Basal.csv

Reactome_2022_high_trans <- read.delim("~/Library/CloudStorage/OneDrive-NorthwesternUniversity/Northwestern/Yi Lab/Manuscript/Code_2ndDraft/Human Skin/Reactome_2022_HIGH_in_trans.txt", quote="")
Reactome_2022_high_trans.sub <- Reactome_2022_high_trans[1:5,]
Reactome_2022_high_trans.sub$Term <- unlist(lapply(Reactome_2022_high_trans.sub$Term , function(x) strsplit(x," R-HSA")[[1]][1]))
Reactome_2022_high_trans.sub$Term <- factor(Reactome_2022_high_trans.sub$Term,levels = rev(Reactome_2022_high_trans.sub$Term))
tiff("~/OneDrive - Northwestern University/Northwestern/Yi Lab/Manuscript/Code_2ndDraft/Human Skin//Reactome_2022_high_trans.tiff", units="in", width=6, height=6, res=300)
ggplot(data=Reactome_2022_high_trans.sub, aes(x=Term, y=-log(Adjusted.P.value))) +
  geom_bar(stat="identity",fill=my.col[2]) + coord_flip()+theme_classic()+theme(axis.title.y = element_blank(),axis.text.y = element_blank(),axis.text = element_text(size = 15,face = "bold"),axis.title  = element_text(size = 15,face = "bold"))
dev.off()


ggplot(data=Reactome_2022_high_trans.sub, aes(x=Term, y=-log(Adjusted.P.value))) +
  geom_bar(stat="identity",fill=my.col[2]) + coord_flip()+theme_classic()+theme(axis.title.y = element_blank(),axis.text = element_text(size = 15,face = "bold"),axis.title  = element_text(size = 15,face = "bold"))+theme_classic()+geom_text(aes(label=Genes), y=c(0),hjust = 0,color="white", size=c(2.5,2.5,2.5,2.5,2.5))+theme(axis.title.y = element_blank())